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[i]  Any  model  of  tides  is  based  on  a  specific  hypothesis  of  how  lagging  depends  on  the 
tidal-flexure  frequency  %.  For  example,  Gerstenkom  (1955),  MacDonald  (1964),  and 
Kaula  (1964)  assumed  constancy  of  the  geometric  lag  angle  8,  while  Singer  (1968)  and 
Mignard  (1979,  1980)  asserted  constancy  of  the  time  lag  At.  Thus  each  of  these  two 
models  was  based  on  a  certain  law  of  scaling  of  the  geometric  lag:  the  Gerstenkom- 
MacDonald-Kaula  theory  implied  that  6  ~  %°,  while  the  Singer-Mignard  theory  postulated 
8  ~  x  ■  The  actual  dependence  of  the  geometric  lag  on  the  frequency  is  more  complicated 
and  is  determined  by  the  rheology  of  the  planet.  Besides,  each  particular  functional 
form  of  this  dependence  will  unambiguously  fix  the  appropriate  form  of  the  frequency 
dependence  of  the  tidal  quality  factor,  Q(x)-  Since  at  present  we  know  the  shape  of  the 
function  Q(x)>  we  can  reverse  our  line  of  reasoning  and  single  out  the  appropriate 
actual  frequency  dependence  of  the  lag,  8(x):  as  within  the  frequency  range  of  our  concern 
Q  ~  x‘\  a  =  0.2-0.4,  then  8  ~  X  °-  This  dependence  turns  out  to  be  different  from  those 
employed  hitherto,  and  it  entails  considerable  alterations  in  the  timescales  of  the  tide¬ 
generated  dynamical  evolution.  Phobos’s  fall  on  Mars  is  an  example  we  consider. 
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1.  Introduction 

[2]  If  a  satellite  is  located  at  a  planetocentric  position  r,  it 
generates  a  tidal  bulge  that  either  advances  or  retards  the 
satellite  motion,  depending  on  the  interrelation  between  the 
planetary  spin  rate  uip  and  the  tangential  part  of  satellite’s 
velocity  v  divided  by  r  =  |r|.  It  is  convenient  to  imagine  (as 
on  Figure  1)  that  the  bulge  emerges  beneath  a  fictitious 
satellite  located  at 

V  =  r  +/,  (1) 

where  the  position  lag  /  is  given  by 

/=  At(up  x  P  -  v).  (2) 

At  is  the  time  lag  between  the  real  and  fictitious  tide¬ 
generating  satellites,  and  the  inclination  and  eccentricity  of 
the  satellite  are  assumed  sufficiently  small. 

[3]  The  fictitious  satellite  is  merely  a  way  of  illustrating 
the  time  lag  between  the  tide-raising  potential  and  the 
distortion  of  the  body.  This  concept  implies  no  new  physics, 
and  is  but  a  convenient  metaphor  employed  to  convey  that 
at  each  instance  of  time  the  dynamical  tide  is  modeled  with 
a  static  tide  where  all  the  time-dependent  variables  are 
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shifted  back  by  At,  i.e.,  (1)  the  moon  is  rotated  back  by 
v  At,  and  (2)  the  attitude  of  the  planet  is  rotated  back  by  u>p 
At.  From  the  viewpoint  of  a  planet-based  observer,  this 
means  that  a  dynamical  response  to  a  satellite  located  at  P  is 
modeled  with  a  static  response  to  a  satellite  located  at  Pf  = 
r  -  At  (v  -  u>p  x  P),  provided  we  take  into  account  only 
the  principal  bulge.  Neglect  of  the  higher  frequencies  is 
permissible  for  small  eccentricities  and  inclination.  In  this 
case  it  is  acceptable  to  use  the  fictitious-satellite  metaphor 
as  a  euphemism  for  bulge  lagging. 

[4]  In  this  paper,  we  intend  to  dwell  on  geophysical 
issues:  the  frequency  dependence  of  the  attenuation  rate 
and  its  consequences.  Hence,  to  avoid  unnecessary  mathe¬ 
matical  complications,  in  the  subsequent  illustrative  example 
we  shall  restrict  ourselves  to  the  simple  case  of  a  tide-raising 
satellite  on  a  near-equatorial  near-circular  orbit.  In  this 
approximation,  the  velocity  of  the  satellite  relative  to  the 
surface  is 

\lZp  x  r  -  \\  =  r\up  -  n\,  (3) 

the  principal  tidal  frequency  is 

X  =  2\ujp-n\,  (4) 

and  the  angular  lag  is 

At  „  At 

6  =  —\upxr~\\=—x,  (5) 

n  being  the  satellite’s  mean  motion,  and  u>p  being  the 
planet’s  spin  rate.  The  factor  of  two  emerges  in  (4)  since 
the  moon  causes  two  elevations  on  the  opposite  sides  of  the 
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Figure  1.  A  planet  and  a  tide-raising  moon.  This  picture 
illustrates  the  case  of  a  satellite  located  below  the 
synchronous  orbit,  so  that  its  mean  motion  n  exceeds 
the  planet’s  spin  rate  u)p,  and  the  _tidal  bulge  is  lagging. 
The  angular  lag  defined  as  6  =  \f\/r  =  y-\^p  x  r  —  v| 
will,  generally,  differ  from  the  absolute  value  of  the  angle 
6 1  subtended  at  the  planet’s  center  between  the  directions 
to  the  satellite  and  the  bulge.  Since  in  our  study  we 
consider  an  example  with  a  small  eccentricity  and 
inclination,  we  make  no  distinction  between  6  and  |<5i|. 


planet.  It  will  also  be  assumed  that  x  At  <C  I ,  for  which 
reason  we  shall  neglect  the  second-order  difference  between 
the  expression  (5)  and  the  angle  subtended  at  the  planet’s 
center  between  the  moon  and  the  tidal  bulge  (rigorously 
speaking,  the  sine  of  the  subtended  angle  is  equal  to  /  x  r/r2). 

[5]  The  starting  point  of  all  tidal  models  is  that  each 
elementary  volume  of  the  planet  is  subject  to  a  tide-raising 
potential,  which  in  general  is  not  periodic  but  can  be 
expanded  into  a  sum  of  periodic  terms.  Within  the  linear 
approximation  introduced  by  Love,  the  tidal  perturbations 
of  the  potential  yield  linear  response  of  the  shape  and  linear 
variations  of  the  stress.  In  extension  of  the  linearity  approx¬ 
imation,  it  is  always  implied  that  the  overall  dissipation 
inside  the  planet  may  be  represented  as  a  sum  of  attenuation 
rates  corresponding  to  each  periodic  disturbance: 

{E)  =  J2  (MXi))  (6) 

i 


where,  at  each  frequency  \u 


(E(X,))  _  _  Epeak(Xi) 
QiXi)  ~  Xi  fi(x'i) 


(V) 


(. . .}  standing  for  averaging  over  flexure  cycle,  E(\j) 
denoting  the  energy  of  deformation  at  the  frequency  ;\y,  and 
Q(xd  being  the  quality  factor  of  the  material  at  this 
frequency.  Introduced  empirically  as  a  means  to  figleaf  our 
lack  of  knowledge  of  the  attenuation  process  in  its  full 
complexity,  the  notion  of  Q  has  proven  to  be  practical  due 
to  its  smooth  and  universal  dependence  upon  the  frequency 
and  temperature.  An  alternative  to  employment  of  the 
empirical  Q  factors  would  be  comprehensive  modeling  of 
dissipation  using  a  solution  of  the  equations  of  motion, 


given  a  rheological  description  of  the  mantle  [. Mitrovica  and 
Peltier,  1992;  Hanyk  et  al.,  1998,  2000;  Moore  and 
Schubert,  2000],  Though  providing  a  valuable  yield  for  a 
geophysicist,  this  comprehensive  approach  may  be  avoided 
in  astronomy,  where  only  the  final  outcome,  the  frequency 
dependence  of  Q,  is  important.  Fortunately,  this  dependence 
is  already  available  from  observations. 

[6]  In  this  paper  we  shall  restrict  ourselves  to  the  simple 
case  of  an  equatorial  or  near-equatorial  satellite  describing  a 
circular  or  near-circular  orbit.  Under  these  circumstances 
only  the  principal  tidal  frequency  (4)  will  matter. 

2.  Quality  Factor  Q  and  the  Geometric 
Lag  Angle  8 

[7]  During  tidal  flexure,  the  energy  attenuation  through 
friction  is,  as  ever,  accompanied  by  a  phase  shift  between 
the  action  and  the  response.  The  tidal  quality  factor  is 
interconnected  with  the  phase  lag  e  and  the  angular  lag  5  via 

Q~l  =  tane  =  tan  2  S  (8) 

or,  for  small  lag  angles, 

Q1  «  e  =  2  S.  (9) 

The  doubling  of  the  lag  is  a  nontrivial  issue.  Many  authors 
erroneously  state  that  Q~  1  is  equal  simply  to  the  tangent  of 
the  lag,  with  the  factor  of  two  omitted.  For  example,  Rainey 
and  Aharonson  [2006]  assume  that  Q~l  is  equal  to  the 
tangent  of  the  geometric  lag.  As  a  result,  they  arrive  at  a 
value  of  Q  that  is  about  twice  larger  than  those  obtained  by 
the  other  teams.  Bills  et  al.  [2005]  used  one  letter,  7,  to 
denote  two  different  angles.  Prior  to  equation  (24)  in  that 
paper,  7  signifies  the  geometric  lag  (in  our  notations,  6)). 
Further,  in  their  equations  (24)  and  (25),  Bills  et  al.  [2005] 
employ  the  notation  7  to  denote  the  phase  lag  (in  our 
notations,  e,  which  happens  to  be  equal  to  2  <5i).  With  this 
crucial  caveat,  Bills  et  aids  [2005]  equation  Q  =  1 /tarry  is 
correct.  The  inaccuracy  in  notations  has  not  prevented  Bills 
et  al.  [2005]  from  arriving  at  a  reasonable  value  of  the 
Martian  quality  factor,  85.58  ±  0.37.  (A  more  recent  study 
by  Lainey  et  al.  [2007]  has  given  a  comparable  value  of 
79.91  ±  0.69.) 

[8]  In  Appendix  A,  section  Al,  we  offer  a  simple 
illustrative  calculation,  which  explains  whence  this  factor 
of  two  stems. 

[9]  Formulae  (8)— (9)  look  reasonable:  the  higher  the 
quality  factor,  the  lower  the  damping  rate  and,  accordingly, 
the  smaller  the  lag.  What  looks  very  far  from  being  OK  are 
the  frequency  dependencies  ensuing  from  the  assertions  of  S 
being  either  constant  or  linear  in  frequency:  the  approach 
taken  by  Gerstenkorn  [1955],  MacDonald  [1964],  and 
Kaula  [1964]  implies  that  Q  ~  y°,  while  the  theory  of 
Singer  [1968]  and  Mignard  [1979,  1980]  yields  Q  ~  %_1, 
neither  option  being  in  agreement  with  the  geophysical  data. 

3.  Dissipation  in  the  Mantle 
3.1.  Generalities 

[10]  Back  in  the  60s  and  70s  of  the  past  century,  when  the 
science  of  low-frequency  seismological  measurements  was 
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yet  under  development,  it  was  widely  thought  that  at  long 
timescales  the  quality  factor  of  the  mantle  is  proportional  to 
the  inverse  of  the  frequency.  This  fallacy  proliferated  into 
planetary  astronomy  where  it  was  received  most  warmly, 
because  the  law  Q~  l/\  turned  out  to  be  the  only  model  for 
which  the  linear  decomposition  of  the  tide  gives  a  set  of 
bulges  displaced  from  the  direction  to  the  satellite  by  the 
same  angle.  Any  other  frequency  dependence  Q{\)  entails 
superposition  of  bulges  corresponding  to  the  separate  fre¬ 
quencies,  each  bulge  being  displaced  by  its  own  angle.  This 
is  the  reason  why  the  scaling  law  Q  ~  l/%,  long  disproved 
and  abandoned  in  geophysics  (at  least,  for  the  frequency 
band  of  our  concern),  still  remains  a  pet  model  in  celestial 
mechanics  of  the  solar  system. 

[11]  Over  the  past  20  a  (where  a  is  years),  considerable 
progress  has  been  achieved  in  the  low-frequency  seismo- 
logical  measurements,  both  in  the  lab  and  in  the  field.  Due 
to  an  impressive  collective  effort  undertaken  by  several 
teams,  it  is  now  a  firmly  established  fact  that  for  frequencies 
down  to  about  ~1  a-1  the  quality  factor  of  the  mantle  is 
proportional  to  the  frequency  to  the  power  of  a  positive 
fraction  a.  This  dependence  holds  for  all  rocks  within  a 
remarkably  broad  band  of  frequencies:  from  several  MHz 
down  to  about  1  a-1. 

[12]  At  timescales  longer  than  1  a,  all  the  way  to  the 
Maxwell  time  (about  100  a),  attenuation  in  the  mantle  is 
defined  by  viscosity,  so  that  the  quality  factor  is,  for  all 
minerals,  well  approximated  with  r\xJM,  where  77  and  M  are 
the  shear  viscosity  and  the  shear  elastic  modulus  of  the 
mineral.  Although  the  values  of  both  the  viscosity  coeffi¬ 
cients  and  elastic  moduli  greatly  vary  for  different  minerals 
and  are  sensitive  to  the  temperature,  the  overall  quality 
factor  of  the  mantle  at  such  long  timescales  still  is  linear  in 
frequency. 

[13]  At  present  there  is  no  consensus  in  the  seismological 
community  in  regard  to  the  time  scales  exceeding  the 
Maxwell  time.  One  viewpoint  (incompatible  with  the  Max¬ 
well  model)  is  that  the  linear  law  Q  ~  x  extends  all  the  way 
down  to  the  zero-frequency  limit  [Karato,  2007],  An 
alternative  point  of  view  (prompted  by  the  Maxwell  model) 
is  that  at  scales  longer  than  the  Maxwell  time  we  return  to 
the  inverse-frequency  law  Q  ~l/%. 

All  in  all,  we  can  write  (using  the  notation  a  for  years) 

For  107  Hz  >  x  >  1  a~*  :  Q  ~  x“, 

with  a  =  0.2  —  0.4  (0.2  for  partial  melts).  (10) 


For  1  a  1  >  x  >  10  2  a  '  :  2  ~  X-  (11) 

For  10~2  a-1  >  \  :  arguably,  it  is  still  Q  ~ 

(Or  maybe  Q  ~  1/x?)  (12) 

Fortunately,  in  practical  calculations  of  tides  in  planets  one 
never  has  to  transcend  the  Maxwell  timescales,  so  the 
controversy  remaining  in  (12)  bears  no  relevance  to  our 
subject.  We  leave  for  a  future  study  the  case  of  synchronous 
satellites,  the  unique  case  of  the  Pluto-Charon  resonance,  or 
the  binary  asteroids  locked  in  the  same  resonance.  Thus  we 


shall  avoid  also  the  frequency  band  addressed  in  (11),  but 
shall  be  interested  solely  in  the  frequency  range  described  in 
(10).  It  is  important  to  emphasize  that  the  positive-power 
scaling  law  (10)  is  well  proven  not  only  for  samples  in  the 
lab  but  also  for  vast  seismological  basins  and,  therefore,  is 
universal.  Hence  this  law  may  be  extended  to  the  tidal 
friction;  validity  of  this  extension  will  be  discussed  in 
section  3.4.1. 

[14]  Below  we  provide  an  extremely  condensed  review  of 
the  published  data  whence  the  scaling  law  (10)  was  derived 
by  the  geophysicists.  The  list  of  sources  will  be  incomplete, 
but  a  full  picture  can  be  obtained  through  the  further 
references  contained  in  the  works  to  be  quoted  below.  For 
a  detailed  treatment,  see  chapter  1 1  of  the  book  by  Karato 
[2007],  which  contains  a  systematic  introduction  into  the 
theory  of  and  experiments  on  attenuation  in  the  mantle. 

3.2.  Circumstantial  Evidence:  Attenuation  in 
Minerals — Laboratory  Measurements 
and  Some  Theory 

[15]  Even  before  the  subtleties  of  solid-state  mechanics 
with  or  without  melt  are  brought  up,  the  positive  sign  of  the 
power  a  in  the  dependence  Q  ~  x  may  be  anticipated  on 
qualitative  physical  grounds.  For  a  damped  oscillator  obey¬ 
ing  z  +  2  /3  z  +  %2  z  =  0,  the  quality  factor  is  equal  to  x/(2/3), 
i.e.,  Q  ~  x- 

[16]  Solid-state  phenomena  causing  attenuation  in  the 
mantle  may  be  divided  into  three  groups:  the  point-defect 
mechanisms,  the  dislocation  mechanisms,  and  the  grain¬ 
boundary  ones. 

[17]  Among  the  point-defect  mechanisms,  most  important 
is  the  transient  diffusional  creep,  i.e.,  plastic  flow  of 
vacancies,  and  therefore  of  atoms,  from  one  grain  boundary 
to  another.  The  flow  is  called  into  being  by  the  fact  that 
vacancies  (as  well  as  the  other  point  defects)  have  different 
energies  at  grain  boundaries  of  different  orientation  relative 
to  the  applied  shear  stress.  This  anelasticity  mechanism  is 
wont  to  obey  the  power  law  Q  ~  x°  with  a  «  0.5. 

[is]  Anelasticity  caused  by  dislocation  mechanisms  is 
governed  by  the  viscosity  law  Q  ~  x  valid  for  sufficiently 
low  frequencies  (or  sufficiently  high  temperatures),  i.e., 
when  the  viscous  motion  of  dislocations  is  not  restrained 
by  the  elastic  restoring  stress.  (At  higher  frequencies  or/and 
lower  temperatures,  the  restoring  force  “pins”  the  defects. 
This  leads  to  the  law  Q  ~  (1  +  t2\2)t~1x~1,  parameter  r 
being  the  relaxation  time  whose  values  considerably  vary 
among  different  mechanisms  belonging  to  this  group.  As 
the  mantle  is  warm  and  viscous,  we  may  ignore  this  caveat.) 

[19]  The  grain-boundary  mechanisms,  too,  are  governed 
by  the  law  Q  ~  though  with  a  lower  exponent:  a  « 
0.2-0. 3.  This  behavior  gradually  changes  to  the  viscous 
mode  (a  =  1)  at  higher  temperatures  and/or  at  lower 
frequencies,  i.e.,  when  the  elastic  restoring  stress  reduces. 

[20]  We  see  that  in  all  cases  the  quality  factor  of  minerals 
should  grow  with  frequency.  Accordingly,  laboratory  meas¬ 
urements  confirm  that,  within  the  geophysically  interesting 
band  of  x,  the  quality  factor  behaves  as  Q  ~  x“  with  a  = 
0.2-0. 4.  Such  measurements  have  been  described  by 
Karato  and  Spetzler  [1990]  and  Karato  [1998],  Similar 
results  were  reported  in  the  works  by  the  team  off  Jackson; 
see,  for  example,  the  paper  by  Tan  et  al.  [1997],  where 
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numerous  earlier  publications  by  that  group  are  also 
mentioned. 

[21]  In  aggregates  with  partial  melt  the  frequency  depen¬ 
dence  of  Q  keeps  the  same  form,  with  a  leaning  to  0.2  [see, 
for  example,  Fontaine  et  al,  2005,  and  references  therein]. 

3.3.  Direct  Evidence:  Attenuation  in  the  Mantle — 
Measurements  on  Seismological  Basins 

[22]  As  we  are  interested  in  the  attenuation  of  tides,  we 
should  be  prepared  to  face  the  possible  existence  of  mech¬ 
anisms  that  may  show  themselves  over  very  large  geolog¬ 
ical  structures  but  not  in  small  samples  explored  in  the  lab. 
No  matter  whether  such  mechanisms  exist  or  not,  we  would 
find  it  safer  to  require  that  the  positive-power  scaling  law 
Q  ~  x'\  even  though  well  proven  in  the  lab,  be  propped 
up  by  direct  seismological  evidence  gathered  over  vast 
zones  of  the  mantle.  Fortunately,  such  data  are  available, 
and  for  the  frequency  range  of  our  interest  these  data 
conform  well  with  the  lab  results.  The  low-frequency 
measurements,  performed  by  different  teams  over  various 
basins  of  the  Earth’s  upper  mantle,  agree  on  the  pivotal  fact: 
the  seismological  quality  factor  scales  as  the  frequency  to 
the  power  of  a  positive  fraction  a;  see,  for  example, 
Mitchell  [1995],  Stachnik  et  al.  [2004],  Shito  et  al.  [2004], 
and  further  references  given  in  these  sources. 

[23]  So  far,  Figure  11  of  Flanagan  and  Wiens  [1998]  is 
the  only  experimental  account  we  know  of,  which  only 
partially  complies  with  the  other  teams’  results.  The  figure 
contains  two  plots  depicting  the  frequency  dependencies  of 
I  IQshaar  and  1  IQcompress-  While  the  behavior  of  both  param¬ 
eters  remains  conventional  down  to  10  1  FIz,  the  shear 
attenuation  surprisingly  goes  down  when  the  frequency 
decreases  to  10“  J  Flz.  Later,  one  of  the  authors  wrote  to 
us  that  “Both  P  and  S  wave  attenuation  becomes  greater  at 
low  frequencies.  The  trend  towards  lower  attenuation  at  the 
lowest  frequencies  in  Figure  11  is  not  well  substantiated” 
(D.  Wiens,  private  communication,  2006).  Hence  the  con¬ 
sensus  on  (10)  stays. 

3.4.  Consequences  for  the  Tides 

3.4.1.  Tidal  Dissipation  Versus  Seismic  Dissipation 

[24]  For  terrestrial  planets,  the  frequency  dependence  of 
the  Q  factor  of  bodily  tides  is  similar  to  the  frequency 
dependence  (10)-(11)  of  the  seismological  Q  factor.  This 
premise  is  based  on  the  fact  that  the  tidal  attenuation  in  the 
mantle  is  taking  place,  much  like  the  seismic  attenuation, 
mainly  due  to  the  mantle’s  nonrigidity.  This  is  a  nontrivial 
fact  because,  in  distinction  from  earthquakes,  the  damping 
of  tides  is  taking  place  both  due  to  nonrigidity  and  self¬ 
gravity  of  the  planet.  Modeling  the  planet  with  a  homoge¬ 
neous  sphere  of  density  p,  rigidity  p,  surface  gravity  g,  and 
radius  R ,  Goldreich  [1963]  managed  to  separate  the  non¬ 
rigidity-caused  and  self-gravity-caused  inputs  into  the  over¬ 
all  tidal  attenuation.  His  expression  for  the  tidal  quality 
factor  has  the  form 


Q  = 


2  gpR\ 

19  p  )’ 


(13) 


Q0  being  the  value  that  the  quality  factor  would  assume 
were  self-gravity  absent.  To  get  an  idea  of  how  significant 


the  self-gravity-produced  input  could  be,  let  us  plug  there 
the  mass  and  radius  of  Mars  and  the  rigidity  of  the  Martian 
mantle.  For  the  Earth’s  mantle,  p  =  65  4-  80  GPa.  Judging 
by  the  absence  of  volcanic  activity  over  the  past  hundred(s) 
of  millions  of  years  of  Mars’  history,  the  temperature  of  the 
Martian  upper  mantle  is  (to  say  the  least)  not  higher  than 
that  of  the  terrestrial  one.  Therefore  we  may  safely 
approximate  the  Martian  p  with  the  upper  limit  for  the 
rigidity  of  the  terrestrial  mantle:  p  =  10 ir  Pa.  All  in  all,  the 
relative  contribution  from  self-gravity  will  look  as 


2_  g  pR 

19  p 


6  7  M2 
767r  p  R4 

2 

j  ^  6.7  x  10  nm3kg  s  2  )  (  6.4  x  1023kg  ) 

40  (  10"  Pa  )  (  3.4  x  10 6m  )* 

5.2  x  10  2, 


(14) 


7  denoting  the  gravity  constant.  This,  very  conservative 
estimate  shows  that  self-gravitation  contributes,  at  most, 
several  percent  into  the  overall  count  of  energy  losses  due  to 
tides.  This  is  the  reason  why  we  extend  to  the  tidal  Q  the 
frequency  dependence  law  measured  for  the  seismic  quality 
factor. 

3.4.2.  Dissipation  in  the  Planet  Versus  Dissipation  in 
the  Satellite 

[25]  A  special  situation  is  tidal  relaxation  toward  the  state 
where  one  body  shows  the  same  side  to  another.  Numerous 
satellites  show,  up  to  librations,  the  same  face  to  their 
primaries.  Among  the  planets,  Pluto  does  this  to  Charon. 
Such  a  complete  locking  is  typical  also  for  binary  asteroids. 
A  gradual  approach  toward  the  synchronous  orbit  involves 
ever-decreasing  frequencies,  eventually  exceeding  the  limits 
of  equation  (11)  and  thus  the  bounds  of  the  present 
discussion.  Mathematically,  this  situation  still  may  be  tack¬ 
led  by  means  of  (10)  until  the  tidal  frequency  x  decreases  to 
I  a  and  then  by  means  of  (11)  while  \  remains  above  the 
inverse  Maxwell  time  of  the  planet’s  material.  Whether  the 
latter  law  can  be  extended  to  longer  timescales  remains  an 
open  issue  of  a  generic  nature  that  is  not  related  to  a  specific 
model  of  tides  or  to  a  particular  frequency  dependence  of  Q. 
The  generic  problem  is  whether  we  at  all  may  use  the 
concept  of  the  quality  factor  beyond  the  Maxwell  time,  or 
should  instead  employ,  beginning  from  some  low  x,  a 
comprehensive  hydrodynamical  model.  In  the  current  work, 
we  address  solely  the  satellite-generated  tides  on  the  planet. 
The  input  from  the  planet-caused  tides  on  the  satellite  will 
be  considered  elsewhere.  The  case  of  Pluto  will  not  be 
studied  here  either.  Nor  shall  we  address  binary  asteroids. 
(Since  at  present  most  asteroids  are  presumed  loosely 
connected,  and  since  we  do  not  expect  the  dependencies 
(10)— (1 1)  to  hold  for  such  aggregates,  our  theory  should 
not,  without  some  alterations,  be  applied  to  such  binaries.) 

[26]  Thus,  since  we  are  talking  only  about  dissipation 
inside  the  planet,  and  are  not  addressing  the  exceptional 
Pluto-Charon  case,  we  may  safely  assume  the  tidal  frequency 
to  always  exceed  1  a  1 .  Thence  (10)  will  render,  for  a  typical 
satellite. 


Q  ~  XQi  whh  a  =  0.2  —  0.4.  (15) 
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Accordingly,  (9)  will  entail 


[31]  Below  we  shall  concentrate  on  the  frequency  depen¬ 
dence  solely. 


8  ~  x  a,  whh  a  =  0.2  —  0.4.  (16) 


[27]  Another  special  situation  is  a  satellite  crossing  a 
synchronous  orbit.  At  the  moment  of  crossing,  the  principal 
tidal  frequency  %  =  2  \ujp  —  n\  vanishes.  As  (9)  and  (10) 
yield  Q  ss  (28)~{  and  Q  ~  x“,  then  we  get  8  ~  \~a  with  a 
positive  a.  Uncritical  employment  of  these  formulae  will 
then  make  one  think  that  at  this  instant  the  lag  8  grows 
infinitely,  a  clearly  nonsensical  result.  The  quandary  is 
resolved  through  the  observation  that  the  bulge  is  lagging 
not  only  in  its  position  but  also  in  its  height,  for  which 
reason  the  dissipation  rate  remains  finite  (M.  Efroimsky, 
Tidal  torques.  A  critical  review  of  some  techniques,  2007; 
available  at  http://arxiv.org/PS_cache/arxiv/pdf/0712/ 
0712.1056vl.pdf).  Since  in  this  paper  we  shall  not  consider 
crossing  of  or  approach  to  synchronous  orbits,  and  since  the 
example  we  aim  at  is  Phobos,  we  shall  not  go  deeper  into 
this  matter  here. 

3.5.  A  Thermodynamical  Aside:  The  Frequency  and 
the  Temperature 

[28]  In  the  beginning  of  the  preceding  subsection  we 
already  mentioned  that  though  the  tidal  0  differs  from  the 
seismic  one,  both  depend  upon  the  frequency  in  the  same 
way,  because  this  dependence  is  determined  by  the  same 
physical  mechanisms.  This  pertains  also  to  the  temperature 
dependence,  which  for  some  fundamental  reason  combines 
into  one  function  with  the  frequency  dependence. 

[29]  As  explained,  from  the  basic  physical  principles,  by 
Karato  [1998,  2007],  the  frequency  and  temperature  depen¬ 
dencies  of  Q  are  inseparably  connected.  The  quality  factor 
can,  despite  its  frequency  dependence,  be  dimensionless 
only  if  it  is  a  function  not  just  of  the  frequency  per  se  but  of 
a  dimensionless  product  of  the  frequency  by  the  typical  time 
of  defect  displacement.  This  time  exponentially  depends  on 
the  activation  energy  A*,  whence  the  resulting  function  is 

Q  ~  [X  / RT)}a .  (17) 

For  most  minerals  of  the  upper  mantle,  A*  lies  within  the 
limits  of  360-540  kJ  mol  .  For  example,  for  dry  olivine  it 
is  about  520  kJ  mol  1 . 

[30]  Thus,  through  fonnulae  (17)  and  (9),  the  cooling  rate 
of  the  planet  plays  a  role  in  the  orbital  evolution  of 
satellites:  the  lower  the  temperature,  the  higher  the  quality 
factor  and,  thereby,  the  smaller  the  lag  8.  For  the  sake  of  a 
crude  estimate,  assume  that  most  of  the  tidal  attenuation  is 
taking  place  in  some  layer,  for  which  an  average  tempera¬ 
ture  T  and  an  average  activation  energy  A*  may  be  intro¬ 
duced.  Then  from  (17)  we  have:  A Q/Q  «  —  a  A*  AT/RT. 
For  a  reasonable  choice  of  values  a  =  0.3  and  A*  =  5.4  x 
105  Jhnol  a  drop  of  the  temperature  from  Ta  =  2000  K  down 
by  AT  =  200  K  will  result  in  A  Q/Q  «  1.  So  a  10% 
decrease  of  the  temperature  can  result  in  an  about  100% 
growth  of  the  quality  factor. 


4.  Formulae 

[32]  The  tidal  potential  perturbation  acting  on  the  tide¬ 
raising  satellite  is 


m)  -  P  (%  •  tf-r}  ?2)  +  p  (m  ■  j  + . . 


(18) 


where  r  =  |r|  and  iy=  \7j\,  while  the  constants  are  given  by 


A2  = 


k2GmR5 

~2 


A  3  = 


k2  GmR1 


(19) 


kn  being  the  Love  numbers.  (For  derivation  of  (18),  see,  for 
example,  MacDonald  [1964,  and  references  therein].) 

[33]  Three  caveats  will  be  appropriate  at  this  point.  First, 
to  (18)  we  should  add  the  potential  due  to  the  tidal  distortion 
of  the  moon  by  the  planet.  That  input  contributes  mainly  to 
the  radial  component  of  the  tidal  force  exerted  on  the  moon, 
and  entails  a  decrease  in  eccentricity  and  semimajor  axis 
[MacDonald,  1964].  Here  we  omit  this  term,  since  our  goal 
is  to  clarify  the  frequency  dependence  of  the  lag.  Second, 
we  acknowledge  that  in  many  realistic  situations  the  k2  and 
sometimes  even  the  k4  term  is  relevant  [Bills  et  al.,  2005]. 
With  intention  to  keep  these  inputs  in  our  subsequent  work, 
here  we  shall  restrict  our  consideration  to  the  leading  term 
only.  Hence  the  ensuing  formula  for  the  tidal  force  will  read 


T=- 


~ik2Gm2R 5 


f 

r 


Jr  f 


O 


if2/*2) 


0(kiGm2R1 /r9 


’ik->Gm1Ri  r _  ,  -  , 

~ - - \rr  -fr  -2r(r-f)[  (20) 

where  we  kept  in  mind  that  0(f2/r2)  =  0{82)  <C  1 .  The  third 
important  caveat  is  that  in  our  further  exploitation  of  this 
formula  we  shall  take  into  account  the  frequency  depen¬ 
dence  of  the  lag  /,  but  not  of  the  parameter  k2.  While  the 
dependence  f(\)  will  be  derived  through  the  interconnec¬ 
tion  of /  with  8(\)  and  therefore  with  Q{\),  the  value  of  k2 
will  be  asserted  constant.  That  the  latter  is  acceptable  can  be 
proven  through  the  following  formula  obtained  by  Darwin 
[1908]  under  the  assumption  of  the  planet  being  a  Maxwell 
body  [see  also  Correia  and  Laskar,  2003]: 


k2  (x)  =  kfluid 


1  +  xW/J 


1  +  (xV//i2)(l  +  19/V(2 gpR)) 


Here  kfluid  is  the  so-called  fluid  Love  number.  This  is  the 
value  that  k2  would  have  assumed  had  the  planet  consisted 
of  a  perfect  fluid  with  the  same  mass  distribution  as  the 
actual  planet.  Notations  p,  p,  g,  and  g  stand  for  the  rigidity, 
mean  density,  surface  gravity,  and  the  radius  of  the  planet. 
For  these  parameters,  we  shall  keep  using  the  estimates 
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from  section  3.4.  The  letter  77  signifies  the  viscosity.  Up  to 
an  order  or  two  of  magnitude,  its  value  may  be 
approximated,  for  a  terrestrial  planet’s  mantle,  with 
ICE2  kg/(m  s).  This  will  yield:  \  ifl/j2  =  (%  1011  s)2, 
wherefrom  we  see  that  in  all  realistic  situations  pertaining  to 
terrestrial  planets  the  frequency  dependence  in  Darwin’s 
formula  will  cancel  out.  Thus  we  shall  neglect  the  frequency 
dependence  of  the  Love  number  k2  (but  shall  at  the  same 
time  take  into  account  the  frequency  dependence  of  Q,  for  it 
will  induce  frequency  dependence  of  all  three  lags). 

[34]  The  interconnection  between  the  position,  time,  and 
angular  lags, 

S  =  y  =  At  '  |  up  X  r  -  v  |  =yX,  (21) 
can  be  equivalently  rewritten  as 


assembled  of  a  homogeneous  medium,  with  a  uniform 
temperature  distribution,  and  if  attenuation  in  this  medium 
were  caused  by  one  particular  physical  mechanism,  then  £ 
would  be  a  relaxation  timescale  associated  with  this 
mechanism  (say,  the  time  of  defect  displacement).  For  a 
realistic  planet,  £  may  be  interpreted  as  a  relaxation  time 
averaged  (in  the  sense  of  Q  =  £"x  ')  over  the  planet’s  layers 
and  over  the  various  damping  mechanisms  acting  within 
these  layers. 

[37]  As  (10)  entails  6  «  1/(2 Q)  =  (1/2)  £~a  \~a,  then 
(21)  necessitates  for  the  position  lag: 

f=r6f=l-r£~aX-ai,  (25) 

and  for  the  time  lag: 

At  =  X-(“+1),  (26) 


„  LAI  - 

f  =  frS  =  r—x  f, 


where 


(22)  £  being  the  planet’s  integral  parameter  introduced  above, 

and  x  being  a  known  function  (4)  of  the  orbital  variables. 
Putting  everything  together,  we  arrive  at 


4  x  r  —  v 
14  x  r-y | 


(23) 


is  a  unit  vector  pointing  in  the  lag  direction. 

[35]  Be  mindful  that  we  assume  the  inclination  and 
eccentricity  to  be  small,  wherefore  the  ratio 


where 


f  = 2 (£xY 


1  —  e2 


1  +  e  cos  v  | up  x  F  —  v| 


1  tr  n-<*  4  x  v  ,  n,  s 

=  2  (£x)  a  1,44;-  +  O(e) 


(27) 


If!  1 

=  At  -  4  x  r  -  v  | 
r  r 

is  simply  the  tangential  angular  lag,  i.e.,  the  geometric  angle 
subtended  at  the  primary’s  center  between  the  moon  and  the 
bulge.  In  the  general  case  of  a  finite  inclination  or/and 
eccentricity,  all  our  formulae  will  remain  in  force,  but  the 
lag  6  will  no  longer  have  the  meaning  of  the  subtended 
angle. 

[36]  At  this  point,  it  would  be  convenient  to  introduce  a 
dimensional  integral  parameter  £  describing  the  overall  tidal 
attenuation  rate  in  the  planet.  The  power  scaling  law 
mentioned  in  section  3  may  be  expressed  as 


Q 


£axa , 


(24) 


where  £a  is  simply  the  dimensional  factor  emerging  in  the 
relation  Q  ~  As  mentioned  in  section  3.5,  cooling  of  the 
planet  should  become  a  part  of  long-term  orbital  calcula¬ 
tions.  It  enters  these  calculations  through  evolution  of  this 
parameter  £.  Under  the  assumption  that  most  of  the  tidal 
dissipation  is  taking  place  in  some  layer,  for  which  an 
average  temperature  T and  an  average  activation  energy  A* 
may  be  defined,  (17)  yields: 


£  =  £ o  exp 


£„  and  Ta  being  the  values  of  the  integral  parameter  and 
temperature  at  some  fiducial  epoch.  The  physical  meaning 
of  the  integral  parameter  £  is  transparent:  if  the  planet  were 


X  =  2  \uip  —  n\.  (28) 

The  time  lag  is,  according  to  (26): 

-(q+1) 

At  =  £  (2  £  | up  -  n|)  .  (29) 

Formulae  (29),  (27),  and  (20)  are  sufficient  to  both  compute 
the  orbit  evolution  and  trace  the  variations  of  the  time  lag. 

5.  Example  of  Phobos’s  Fall  to  Mars 

[38]  As  an  illustrative  example,  let  us  consider  how  the 
realistic  dependence  Q(x)  alters  the  life  time  left  for 
Phobos.  We  shall  neglect  the  fact  that  Phobos  is  close  to 
its  Roche  limit,  and  may  be  destroyed  by  tides  prior  to  its 
fall.  We  also  shall  restrict  the  dynamical  interactions  be¬ 
tween  Phobos  and  Mars  to  a  two-body  problem  disturbed 
solely  with  the  tides  raised  by  Phobos  on  Mars.  Thus  we 
shall  omit  all  the  other  perturbations,  like  the  Martian  non¬ 
sphericity  and  precession,  or  the  pull  exerted  upon  Phobos 
by  the  Sun,  the  planets,  and  Deimos.  If,  along  with  these 
simplifications,  we  assume  the  eccentricity  and  inclination 
to  be  small,  then  we  shall  be  able  to  describe  the  evolution 
of  the  semimajor  axis  by  means  of  the  following  equation 
[Kaula,  1964,  p.  677,  formula  (41)]: 

da  _ _ 3  k2  Rs  Gm 

dt  Q  4 G(M0  +  m)  a11/2 

with  Ma  and  m  denoting  the  masses  of  Mars  and  Phobos. 
This  equation  can  be  solved  analytically,  provided  the 
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Time  (years) 


Figure  2.  Evolution  of  Phobos’s  semimajor  axis,  as  predicted  by  different  models.  The  lines  (from  left 
to  right)  correspond  to  the  Singer-Mignard  model,  to  the  Gerstenkom-McDonald-Kaula  model,  and  to  the 
realistic  rheology  with  a  =  0.3,  with  a  =  0.2,  and  with  a  =  0.4. 


quality  factor  Q  is  set  constant  (as  done  by  Kaula  [1964]). 
The  solution  is: 


a{t) 


(_  39h#Gm 
V  2  Q^G{M  +  m) 


(31) 


aQ  =  a(t)|r=o  being  the  initial  value  of  Phobos’s  semimajor 
axis. 

[39]  Unfortunately,  neither  our  model  (wherein  Q  is  given 
by  (24))  nor  the  Singer-Mignard  model  (with  Q  scaling  as 
1  /%)  admit  such  an  easy  analytical  solution.  This  compels 
us  to  rely  on  numerics.  In  a  (quasi)inertial  frame  centered 
at  Mars,  the  equation  of  motion  looks: 


<£7 

dt2 


G(M0  +  m)7  3k2G(M0  +  m)mRs  r  -  2  n  /  -at 

?  >M0  V  -V  \ 

(32) 


Naturally,  its  right-hand  side  consists  of  the  principal,  two- 
body  contribution  -  G(Ma  +  in)7l£  and  the  disturbing  tidal 
force  given  by  (20).  It  should  be  noted,  however,  that  we  did 
not  bring  in  here  all  the  terms  from  (20).  Following  Mignard 
[1980],  we  retain  in  (32)  only  the  perturbing  terms 
dependent  on  /.  The  other  perturbing  term  (the  first  term 
on  the  right-hand  side  of  (20))  is  missing  in  (32),  because  it 
provides  no  secular  input  into  the  semimajor  axis’s  evolution. 
For  a  proof  of  this  fact  see  Appendix  A,  section  A3. 

[40]  Phobos’s  orbital  motion  may  be  described  with  the 
planetary  equations  in  the  Euler-Gauss  form.  In  assumption 


of  i  and  e  being  small,  the  problem  conveniently  reduces  to 
one  equation: 


da  _  2  S 
dt  n  ’ 


(33) 


where  5  is  the  tidal  acceleration  given  by  the  second  term 
of  (32)  projected  onto  the  instantaneous  direction  of  the 
satellite  motion.  This  direction  is  defined  by  the  unit 
vector  (H  x  7)/\H  x  7  |,  where  H  denotes  the  angular- 
momentum  vector.  Thence  the  said  projection  reads: 


S  = 


3k2G(Mc  +  m)?nR5 


r'°M„ 


[— r1  IS.t(u:p  x  7  —  v)  +  27 At (7  ■  v)] 


(H  x  7) 
I H  x  7\ 


(34) 


Mignard  [1980]  gives  the  appropriate  expression  with  a 
wrong  sign.  This  is  likely  to  be  a  misprint,  because  the 
subsequent  formulae  in  his  paper  are  correct. 

J41]  In  assumption  of  i  and  e  being  negligibly  small, 
(H  x  7)  ■  (c dp  x  7)  can  be  approximated  with  n  a4  up, 
whereafter  (34)  gets  simplified  to 

3k2R5G(M0  +  m)mAt  ,  , 

5= - —7 - (»-<*),  05) 


substitution  whereof  into  (33)  entails  the  following  equation 
to  integrate: 


da 

dt 


6k2R5nmAt 

Mccd 


(«-Wp). 


(36) 
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Figure  3.  Evolution  of  At  computed  for  three  different  models.  The  upper  horizontal  line  corresponds 
to  the  Singer-Mignard  model  wherein  At  is  set  to  be  constant.  The  other  four  lines  correspond  to  the 
Gerstenkom-MacDonald-Kaula  model  and  to  the  realistic  rheology  with  a  =  0.3,  with  a  =  0.2,  and  with 
a  =  0.4. 


Our  computational  scheme  was  based  on  the  numerical 
integrator  RA15  offered  by  Everhart  [1985],  The  initial 
value  of  a,  as  well  as  the  values  of  all  the  other  physical 
parameters  entering  (36),  were  borrowed  from  Lainey  et  al. 
[2007].  These  included  an  estimate  of  0.6644  min  for  the 
present-day  time  lag  At. 

[42]  Four  numerical  simulations  were  carried  out.  One  of 
these  implemented  the  Singer-Mignard  model  with  a  tidal- 
frequency-independent  At.  The  other  three  integrations 
were  performed  for  the  realistic  frequency  dependence 
(29),  with  a  =  0.2,  0.3  and  0.4. 

[43]  To  find  the  integral  parameter  £  emerging  in  (26),  we 
used  the  present-day  values  of  At  and  n.  The  resulting  values 
of£  were  found  to  be  1201  x  105  day  rad”1, 81028  day  rad”1 
and  2104  day  rad-  'for  a  =  0.2,  0.3  and  0.4,  respectively. 
We  did  not  take  into  account  the  strong  temperature  depen¬ 
dence  of  £,  leaving  this  interesting  topic  for  discussion 
elsewhere. 

[44]  Simultaneous  numerical  integration  of  equations  (29) 
and  (36)  results  in  plots  presented  on  Figures  2  and  3.  The 
first  of  these  pictures  shows  the  evolution  of  Phobos’s 
semimajor  axis  from  its  present  value  until  the  satellite 
crashes  on  Mars,  having  descended  about  6000  km.  The 
leftmost  curve  reproduces  the  known  result  that,  according 
to  the  Singer-Mignard  model  with  a  constant  At,  Phobos 
should  fall  on  Mars  in  about  29  Ma.  The  next  curve  was 
obtained  not  numerically  but  analytically.  It  depicts  the 
analytical  solution  (31)  available  for  the  Gerstenkom- 
MacDonald-Kaula  model  with  a  constant  Q,  and  demon¬ 
strates  that  this  model  promises  to  Phobos  a  longer  age, 
38  Ma.  The  three  curves  on  the  right  were  obtained  by 
numerical  integration  of  (29)  and  (36).  They  correspond  to 


the  realistic  rheology  with  a  equal  to  0.2,  0.3  and  0.4.  It  can 
be  seen  that  within  the  realistic  model  Phobos  is  expected  to 
survive  for  about  40-43  Ma,  dependent  upon  the  actual 
value  of  a  of  the  Martian  mantle.  This  is  about  15  Ma 
longer  than  within  the  Singer-Mignard  model  widely  ac¬ 
cepted  hitherto. 

[45]  (In  the  paragraph  after  his  formula  (18),  Mignard 
[1981]  states  that  “Phobos  will  end  its  life  in  about  36 
million  years.”  Mignard  arrived  at  that  number  by  using  an 
old  estimate  of  20  deg/ca2  for  the  initial  tidal  acceleration. 
Later  studies,  for  example,  Jacobson  et  al.  [1989]  and 
Lainey  et  al.  [2007],  have  shown  that  this  value  should  be 
increased  to  25.4  deg/ca2.  It  is  for  this  reason  that  our 
simulation  based  on  the  Singer-Mignard  model  gives  not  36 
but  only  29  Ma  for  Phobos’s  remaining  lifetime.) 

[46]  The  difference  between  the  three  scenarios  shown  on 
Figure  2  stems  from  the  different  rate  of  evolution  of  the  lag 
At  in  the  three  theories  addressed.  Within  the  Singer- 
Mignard  formalism,  At  stays  unchanged  through  the  de¬ 
scent.  As  can  be  seen  from  formula  (29),  this  is  equivalent 
to  setting  a  =  —  1 ,  an  assertion  not  supported  by  geophys¬ 
ical  data.  Within  the  Gerstenkom-MacDonald-Kaula  model, 
the  time  lag  is  subject  to  a  gradual  decrease  described  by  the 
formula 

At  =  arCten(1/g)  (37) 

2  | n  -  ujp\ 

under  the  assumption  that  Q  is  constant  and  is  equal  to  its 
present-day  value  Q  =  79.91  determined  by  Lainey  et  al. 
[2007].  Comparison  of  (37)  with  (29)  reminds  us  of  the  simple 
fact  that,  in  terms  of  our  model,  Gerstenkom-MacDonald- 
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Kaula’s  theory  corresponds  to  the  choice  of  a  =  0,  a  choice 
which  is  closer  to  the  realistic  rheology  than  the  Singer- 
Mignard  model. 

[47]  In  the  realistic  model,  a  is  positive  and  assumes  a 
value  of  about  0.2-0. 4.  As  a  result,  the  time  lag  is  gradually 
decreasing.  However  this  decrease  looks  different  from  that 
in  Gerstenkom-MacDonald-Kaula’s  model;  for  their  com¬ 
parison,  see  Figure  3. 

6.  Conclusions 

[48]  As  the  tidal  angular  lag  6  is  inversely  proportional  to 
the  tidal  Q  factor,  the  actual  frequency  dependence  of  both  S 
and  At  is  unambiguously  defined  by  the  frequency  depen¬ 
dence  of  Q.  While  in  the  Gerstenkom-MacDonald-Kaula 
theory  of  tides  the  geometric  lag  is  assumed  frequency- 
independent,  in  the  Singer-Mignard  theory  it  is  the  time  lag 
that  is  spared  of  frequency  dependence.  However,  neither  of 
these  two  choices  conform  to  the  geophysical  data. 

[49]  We  introduce  a  realistic  tidal  model,  which  permits 
the  quality  factor  and  therefore  both  the  angular  lag  6  and 
the  time  lag  At  to  depend  on  the  tidal  frequency  The 
quality  factor  is  wont,  according  to  numerous  studies,  to 
obey  the  law  Q  ~  \a,  where  a  lies  within  0.2-0. 4.  This 
makes  the  time  lag  At  not  a  constant  but  a  function  (26)  of 
the  principal  tidal  frequency  and,  through  (29),  of  the  orbital 
elements  of  the  satellite.  The  same  pertains  to  the  angular 
lag  6. 

[50]  Using  these  tidal-frequency  dependencies  for  the 
time  and  angular  lags,  along  with  the  recently  updated 
values  of  the  Martian  parameters,  we  explored  the  future 
of  Phobos,  taking  into  account  only  the  tides  raised  by 
Phobos  on  Mars,  but  not  those  caused  by  Mars  on  Phobos. 
Our  integration  shows  that  Phobos  will  fall  on  Mars  in  40- 
43  Ma  from  now.  It  is  up  to  50%  longer  than  the  estimate 
stemming  from  the  Singer-Mignard  model  employed  in  the 
past.  This  demonstrates  that  the  currently  accepted  time- 
scales  of  dynamical  evolution,  deduced  from  old  tidal 
models,  should  be  reexamined  using  the  actual  frequency 
dependence  of  the  lags. 

Appendix  A 

[51]  The  goal  of  this  appendix  is  threefold.  First,  we 
remind  the  reader  why  in  the  first  approximation  the  quality 
factor  is  inversely  proportional  to  the  phase  lag.  Second,  we 
explain  why  the  phase  lag  is  twice  the  geometric  lag  angle, 
as  in  formulae  (8)— (9)  above.  While  a  comprehensive 
mathematical  derivation  of  this  fact  can  be  found  elsewhere 
(see  the  unnumbered  formula  between  equations  (29)  and 
(30)  on  p.  673  of  Kaula  [1964]),  here  we  illustrate  this 
counterintuitive  result  by  using  the  simplest  setting.  Third, 
we  justify  our  neglect  of  the  first  term  in  (20). 

Al.  The  Case  of  a  Near-Circular  Near-Equatorial 
Orbit 

[52]  Consider  the  simple  case  of  an  equatorial  moon  on  a 
circular  orbit.  At  each  point  of  the  planet,  the  tidal  potential 
produced  by  this  moon  will  read 

W  =  Wa  cosyf,  (Al) 


the  tidal  frequency  being  given  by 

X  =  2\n  —  iop\.  (A2) 

Let  g  denote  the  free-fall  acceleration.  An  element  of  the 
planet’s  volume  lying  beneath  the  satellite’s  trajectory  will 
then  experience  a  vertical  elevation  of 

W 

£  =  — -  cos(xl  —  26).  (A3) 

g 

Accordingly,  the  vertical  velocity  of  this  element  of  the 
planet’s  volume  will  amount  to 

u  =  C=-x—  sin(xt  -  28) 
g 
ffr 

=  —  \  — -  (sin  xt  cos  28  —  cosxlsin2i5).  (A4) 

The  expression  for  the  velocity  has  such  a  simple  form 
because  in  this  case  the  instantaneous  frequency  \  is 
constant.  The  satellite  generates  two  bulges  (on  the  facing 
and  opposite  sides  of  the  planet),  so  each  point  of  the 
surface  is  uplifted  twice  through  a  cycle.  This  entails  the 
factor  of  two  in  the  expressions  (A2)  for  the  frequency.  The 
phase  in  (A3),  too,  is  doubled,  though  the  necessity  of  this  is 
less  evident.  Let  x  signify  a  position  along  the  equatorial 
circumference  of  the  planet.  In  the  absence  of  lag,  the  radial 
elevation  at  a  point  x  would  be 


g 

v  being  the  velocity  of  the  satellite’s  projection  on  the 
ground,  R  being  the  planet’s  radius,  and  a  being  simply  | n  — 
uip\  because  we  are  dealing  with  a  circular  equatorial  orbit. 
The  value  of  k  must  satisfy 

kv  =  2s,  i.e.,  kv  =  x, 

to  make  sure  that  at  each  x  the  ground  elevates  twice  per  an 
orbital  cycle.  The  above  two  formulae  yield: 

kR  =  2. 

In  the  presence  of  lag,  all  above  stays  in  force,  except  that 
the  fonnula  for  radial  elevation  will  read: 

W 

C= — -  cosk(x  —  vt  +  D),  where  D  =  R8, 
g 

D  being  the  linear  lag,  and  6  being  the  angular  one.  Since 
k  v  =  2,  we  get: 

cos[£(x  —  vt  +  Rd, )]  =  cos[foc  —  kv  t  +  kRd] 

=  cosjfct  —  (kvt  —  2d)], 
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so  that,  at  some  fixed  point  (say,  at  x  =  0)  the  elevation 
becomes: 


C(f)  =  — -  cos(kvt  —  2d), 
g 


We  see  that,  while  the  geometric  lag  is  6,  the  phase  lag  is 
double  thereof. 

[53]  The  energy  dissipated  over  a  time  cycle  T=  2tt/\,  per 
unit  mass,  will,  in  neglect  of  horizontal  displacements,  be 


AE„,e  = 


f ■ 


dW 

l)r 


dt=-\-x 


W0\  dw„ 

dr 


rt=T 


■  /  cos  \'f(sin  \t  cos  26  —  cos  \t  sin  26)  dt 

Jt= 0 

w„  dw„  .  .  1  rxt=2* 


=  -X 


g  dr 
W„  dWa 


-  sin  26 


d 


X  J\t= 0 


cos  \t  d(\t) 


g  dr 


7rsin2i5, 


(A5) 


while  the  peak  energy  stored  in  the  system  during  the  cycle 
will  read: 


E peak 


I' 


r/4  /  dw 


dr 


W0\  dWQ 


dr 


/  cosx^(sinx^cos2^  —  cos  sin  2^)  dt 

Jt= o 

cos  26  f  ^=7r/^ 

cos  \t  sin \t  d(\t) 


=  2(7- 


W„  dWa 
dr 


X  J  \t- 


r-\t=n 

Jyt= 0 


12  6  rx‘=n/ 2 


rx'=  . 

/  cos2  \t  d(xt) 

J\t= o 


X  J\t= ( 


Wo  dWq 

g  dr 


cos  2  6 - sin  26 


(A6) 


whence 


Q- 1  =  ~AE cycle  _  1  7T  sin  26 


2n Epeat  2tt  j  cos  26  —  |  sin  26 


:  tan  26.  (A7) 


The  above  formulae  were  written  down  in  neglect  of 
horizontal  displacements,  approximation  justified  below  in 
the  language  of  continuum  mechanics. 

A2.  On  the  Validity  of  our  Neglect  of  the 
Horizontal  Displacements 

[54]  In  our  above  derivation  of  the  interrelation  between 
Q  and  6,  we  greatly  simplified  the  situation,  taking  into 
account  only  the  vertical  displacement  of  the  planetary 
surface,  in  response  to  the  satellite’s  pull.  Here  we  shall 
demonstrate  that  this  approximation  is  legitimate,  at  least  in 
the  case  when  the  planet  is  modeled  with  an  incompressible 
and  homogeneous  medium. 

[55]  As  a  starting  point,  recall  that  the  tidal  attenuation 
rate  within  a  tidally  distorted  planet  is  well  approximated 
with  the  work  performed  on  it  by  the  satellite.  (A  small 
share  of  this  work  is  being  spent  for  decelerating  the  planet 
rotation.)  Thus  we  shall  write  down  the  attenuation  rate  as 


E 


J  pVVWd3x, 


(A8) 


where  p,  V  and  W  are  the  density,  velocity,  and  tidal 
potential  inside  the  planet.  To  simplify  this  expression,  we 
shall  employ  the  equality 

pVVW  =  V  ■  (pVW)  —  WV  ■  V p  —  WV  ■  (pV) 

=  V  ■  (pVVW)  —  WV  ■  Vp+  W-^.  (A9) 

For  a  homogeneous  and  incompressible  primary,  both  the  v 
W  V  p  and  dp/d  t  terms  are  nil,  wherefrom 

E  =  -  j  pWV  ■  nd3x,  (A10) 

n  being  the  outward  normal  to  the  surface  of  the  planet.  We 
immediately  see  that,  within  the  hydrodynamical  model,  it 
is  only  the  radial  elevation  rate  that  matters. 

[56]  Now  write  the  potential  as  W=  W0  cos(x  t).  Since  the 
response  is  delayed  by  At,  the  surface-inequality  rate  will 
evolve  as  V  ■  n  ~  sin[;\;  (/  —  At)].  All  the  rest  will  then  be 
as  in  section  Al. 

A3.  On  the  Validity  of  our  Neglect  of  the 
Nondissipative  Tidal  Potential 

[57]  The  right-hand  side  of  equation  (32)  consists  of  the 
principal  part,  —  G  ( Ma  +  m)  r/r3,  and  tidal  perturbation 
terms.  These  are  the  second  and  third  terms  from  the  right- 
hand  side  of  (20),  terms  that  bear  a  dependence  on  /  and, 
therefore,  on  At.  The  first  term  from  the  right-hand  side  of 
(20)  lacks  such  a  dependence  and  therefore  is  omitted  in 
(32).  The  term  was  dropped  because  it  would  provide  no 
secular  input  into  the  history  of  the  semimajor  axis.  Here  we 
shall  provide  a  proof  of  this  statement. 

[58]  The  omitted  term  corresponds  to  a  potential 
[Mignard,  1980]: 


U0 


ki{M0  +  m)GmR 5 
2Mar6 


(All) 


From  the  physical  standpoint,  Uq  models  the  effect  of  the 
tidal  bulges,  assuming  their  direction  to  coincide  with  that 
toward  the  tide-raising  satellite.  This  potential  entails  no 
angular-momentum  exchange,  and  therefore  yields  no 
secular  effect  on  the  semimajor  axis.  To  prove  this,  let  us 
decompose  this  potential  into  a  series  over  the  powers  of  e. 
This  will  require  of  us  to  derive  the  expression  of  (a/r)6. 
Starting  out  with  the  well  known  development 


-  =  1  +  2  Jp{pe)  cos )pM) 

r  p= 1 


1 


=  1  +  ecos(M)  +  e~  cos(2M)  +  cos(3M)  —  -  cos )M)je3 

+  -  cos(2M)  +  -  cos(4M)^  eA 

(  81  .  .  1  ,  ,  625  ,  A  , 

+  V  128  C0S(3M)  +  192  c°dM)  +  3g4  cos(5 M)\  e 

81  1  ,  ,  16  ,  A 

—  cos(6M)  +  —  cos(2M)  —  —  cos(4M)  Jeb  +  ... , 

(A12) 
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one  can 

\  6 


0‘ 


-- 1 


(y+ycos(2M))e2 


ie3 


arrive  al 

6ecos(M) 

/ 107  ,  ,  117  .  A 

-  (  cos(3M)  +  cos(M)  j 

/ 197  ,  x  101  ,  105 

-  I  cos(4M)  +  ~y~  cos(2 M)  + 

/5 157  5529  2721 

-  cos(3M)  +  —  cos(5M)  +  -L-  cos (M) 

/  4839  ,  x  N 

-  (  cos(4M)  +  129  cos(2M) 

_  525  +  732  cqs(6M)\  6+  ^  ( 


whose  average  over  the  mean  anomaly  looks  like: 


Hence  the  averaged  potential  will  become: 


(A14) 


(Uo) 


k2{M0  +  m)GmRs  (  15  2  105  4 

- ; -  14 - e~  H - e 

2  M0a6  V  2  4 


(A15) 


In  the  case  of  Phobos,  the  terms  of  order  0(e2)  may,  in  the 
first  approximation,  be  neglected.  This  means  that  out  of  the 
six  Lagrange-type  planetary  equations  the  first  five  will,  in 
the  first  order  of  e,  stay  unperturbed,  and  therefore  the 
elements  a,  e,  ui,  i,  O  will,  in  the  first  order  over  e,  remain 
unchanged.  The  Lagrange  equation  for  the  longitude  will  be 
the  only  one  influenced  by  (70-  That  equation  will  assume 
the  form: 


dL  _  2  d(U0) 

dt  na  da 


which  gives 


L  =  nt  + 


6k2{M0  +  m)GmR 5 
M„na 8 


(A16) 


(A17) 


We  see  that  in  the  first  order  of  e  the  only  secular  effect 
stemming  from  the  potential  Uq  is  a  linear  in  time  evolution 
of  the  longitude. 
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